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16.  Abstract 

The  Coast  Guard  uses  ocean  surface  velocity  fields  as  input  to  trajectory  models  in  search  planning  to 
determine  search  area  movement.  Statistical  techniques,  called  objective  analysis,  are  used  by 
meteorologists  and  oceanographers  to  combine  observed  data,  which  are  irregularly  distributed  in  space  and 
time,  to  make  conclusions  about  the  total  variation  of  the  variables  within  a  region  of  interest.  The  purpose  of 
the  Phase  I  research  was  to  study  the  feasibility  of  adapting  objective  analysis  techniques  for  use  with  Coast 
Guard  data.  The  following  research  was  completed.  Objective  analysis  techniques  were  used  to  combine 
irregularly  spaced  and  temporally  observed  ocean  current  velocity  data  in  order  to  determine  an  optimally 
estimated  velocity  field  and  associated  error  field.  The  technique  was  applied  to  a  drifting  buoy  data  set 
collected  by  the  Coast  Guard.  A  numerical  experiment  was  also  conducted  to  examine  the  number  of  buoys 
required  to  determine  the  velocity  field  in  a  typical  Search  and  Rescue  domain.  Techniques  for  combining  two 
velocity  fields  were  also  investigated. 

The  research  findings  were  that  objective  analysis  techniques  can  be  used  to  provide  the  Coast  Guard  with 
accurate  estimates,  in  real  time,  of  ocean  surface  velocity  and  associated  error. 

Potential  applications  of  this  research  are  in  support  of  Coast  Guard  Search  and  Rescue,  Ice  Patrol  and 
Environmental  Protection  missions. 
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1.0  INTRODUCTION 


The  Coast  Guard  uses  ocean  surface  velocity  fields  from  various  ocean  circulation  models  to  pro¬ 
vide  input  into  trajectory  models  which  are  used  in  search  planning  to  determine  search  area 
movement.  Numerical  techniques,  called  objective  analysis,  are  used  by  oceanographers  to  com¬ 
bine  observed  data,  irregularly  distributed  in  space  and  time,  to  determine  the  total  variation  of  the 
variables  in  the  region  of  interest. 

The  purpose  of  this  project  is  to  study  the  feasibihty  of  adapting  objective  analysis  techniques  for 
use  with  Coast  Guard  observation  data  in  order  to  provide  a  significant  improvement  in  Coast 
Guard  search  effectiveness  in  a  variety  of  mission  areas. 

The  following  tasks  were  completed: 

Use  objective  analysis  techniques  to  combine  irregularly  spaced  and  tem¬ 
porally  observed  ocean  current  velocity  data  in  order  to  determine  an  opti¬ 
mally  estimated  velocity  field,  and 

Develop  techniques  to  combine  the  regularly  gridded  velocity  data  output 
from  two  or  more  ocean  circulation  models  to  obtain  a  single  optimally 
estimated  velocity  field. 

In  this  report,  we  present  the  results  of  the  Phase  I  effort. 

Organization  of  Report.  The  report  is  organized  in  the  following  manner. 

Coast  Guard  applications  of  ocean  surface  velocity  fields  are  discussed  in  Section  2.  A  descrip¬ 
tion  of  the  mathematical  theory  underlying  Objective  Analysis  (OA)  is  given  in  Section  3.  A 
computer  implementation  of  the  OA  algorithm  to  estimate  ocean  surface  velocity  fields  is 
described  in  Section  4. 

The  Dataset  used  to  illustrate  the  apphcation  of  Objective  Analysis  is  described  in  Section  5. 
Examples  of  an  Objective  Analysis  derived  velocity  field  are  given  in  Section  6.  A  method  to 
combine  velocity  fields  from  two  or  more  models  is  given  in  Section  7.  An  illustration  of  the  pro¬ 
jection  of  an  area  of  uncertainty  using  an  OA  velocity  field  is  given  in  Section  8.  Conclusions  are 
presented  in  Section  9. 
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2.0  COAST  GUARD  APPLICATIONS  OF  OCEAN  SURFACE 
VELOCITY  FIELDS 


The  Coast  Guard  is  the  primary  federal  agency  with  maritime  authority  for  the  United  States.  The 
Coast  Guard’s  four  main  missions  are  maritime  law  enforcement,  maritime  safety,  environmental 
protection  and  national  security. 

Specific  Coast  Guard  missions  that  will  benefit  from  this  research  are  Search  and  Rescue,  Ice 
Patrol,  and  Environmental  Protection.  Central  to  the  Coast  Guard’s  research  efforts  in  support  of 
these  missions  is  the  need  to  model  accurately  the  ocean  current  in  order  to  calculate  the  move¬ 
ment  of  search  areas. 

Search  and  Rescue  (SAR).  The  Coast  Guard  operates  the  world’s  largest  SAR  organization 
consisting  of  a  nationwide  system  of  boats,  aircraft,  cutters,  and  rescue  coordination  centers  on 
24-hour  alert  ready  to  rescue  persons  in  peril  at  sea.  The  objective  of  the  Coast  Guard  SAR  pro¬ 
gram  is  to  minimize  loss  of  life,  personal  injury  and  property  damage  on  the  high  seas  and  in  all 
U.S.  waters. 

In  1992,  a  total  of  52,645  SAR  cases  were  conducted  with  5,547  lives  saved  and  $538M  in  prop¬ 
erty  loss  prevented. 

Ice  Patrol.  The  Coast  Guard  conducts  the  International  Ice  Patrol  to  observe  and  chart  the  posi¬ 
tion  and  movements  of  icebergs. 

Environmental  Protection.  The  Coast  Guard  is  the  lead  agency  to  respond  to  the  threat  of  pollu¬ 
tion  in  the  coastal  zone.  Their  mission  is  to  minimize  damage  caused  by  spills  of  oil  or  hazardous 
substances  and  to  assist  in  pollution  response  planning  efforts.  In  1992,  the  Coast  Guard 
responded  to  8,700  oil  pollution  incidents. 
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3.0  OBJECTIVE  ANALYSIS 


3.1  INTRODUCTION 

Objective  Analysis  is  a  numerical  procedure  for  combining  a  number  of  atmospheric  or  oceano¬ 
graphic  observations,  which  may  be  irregularly  distributed  in  space  and  time,  in  order  to  draw 
conclusions  about  the  total  variation  of  meteorological  or  oceanographic  variables  in  a  specific 
region  (Reference  [a]). 

Numerical  procedures  that  have  been  used  for  objective  analysis  include  polynomial  interpola¬ 
tion,  statistical  interpolation,  successive  correction  methods,  variational  methods,  and  spectral 
analysis  methods.  A  study  conducted  by  Gandin,  discussed  in  Reference  [a],  indicated  that  when 
the  data  available  are  dense,  these  techniques  all  give  similar  results.  When  the  data  are  sparse, 
the  statistical  interpolation  method  provided  the  best  performance.  The  objective  analysis  tech¬ 
nique  used  in  this  report  is  a  statistical  interpolation  procedure  which  provides  a  minimum  vari¬ 
ance  estimate. 

The  description  of  objective  analysis  in  this  section  is  taken  from  References  [a]  and  [b].  The 
mathematical  theory  for  objective  analysis  of  a  scalar  field  (i.e.,  a  single  variable  specified  at  each 
space-time  location)  is  given  in  Section  3.2.1.  These  results  are  extended  in  Section  3.2.2  to  apply 
to  a  vector  field  with  several  variables  specified  at  each  space-time  point,  and  where  different  field 
variables  may  be  cross-correlated. 


3.2  MATHEMATICAL  THEORY 

3.2.1  Scalar  Objective  Analysis 

To  explain  the  theory  behind  objective  analysis,  we  introduce  the  following  notation  and  statisti¬ 
cal  concepts. 

The  quantity  under  observation  vrill  be  referred  to  as  “the  field”,  and  denoted  by  the  symbol/.  In 
addition,  let 

denote  the  true  field  value  at  position  i  and  time  t,  (la) 

the  observed  field  value  at  position  i  and  time  t,  and  (lb) 

gj.,,  the  observational  error  at  position  i  and  time  t.  (Ic) 

These  quantities  are  related  by  the  equation 

(2) 
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For  a  continuous  variable,  a  mean  value  over  time  will  be  denoted 


T 

For  a  discrete  variable,  the  mean  value  over  n  observations  will  be  denoted 

l4,- 

(3b) 

The  deviation  of  a  field  from  its  mean  (sometimes  referred  to  as  an  anomaly)  will  be  represented 
by 


/'<,  =4-V„>' 


(4) 


The  mean  value  over  time  of  the  square  of  the  deviation  from  the  mean  is  known  as  the  variance, 

.  Its  square  root  is  called  the  standard  deviation,  a . .  The  covariance  between  the  field  at  posi¬ 
tion  i  and  position  j  is  denoted  by 

Because  the  true  field  value  is  usually  not  known,  the  concept  of  the  anomaly  must  be  generalized 
by  allowing  it  to  be  the  deviation  of  an  observation  from  a  preliminary  estimate  of  the  true  field 
value.  A  convenient  preliminary  estimate  is  the  weighted  mean  value  of  temporally  and  spatially 

nearby  points.  The  symbol  fP^  will  denote  the  preliminary  field  value  which  is  given  by 


fobs 

Jjt 


j=i 


(6) 


The  weights,  l-j,  in  Equation  (6)  typically  depend  on  the  spatial  and  temporal  distance  of  the 
observations  from  the  point  of  interest  and  may  include  a  dependence  on  the  amount  of  measure¬ 
ment  error  associated  with  an  observation. 
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In  the  following,  the  field  anomalies  and  covariances  will  be  taken  with  respect  to  the  preliminary 
estimate  of  the  field,  rather  than  its  mean  value  over  time.  Thus,  Equations  (4)  and  (5)  are 
replaced  by 


fu=ftrfu’^  P) 

"•tj  -  Wu-fS) 

=  Vu'fiP-  ® 

It  is  assumed  that  /;/  has  mean  0,  that  is,  (fi/)  =  0. 

After  the  preliminary  field  value  is  determined,  a  minimum  variance  estimate  of  the  field  at  a 
point  can  be  calculated.  A  minimum  variance  estimate  is  one  whose  value  will  on  average  deviate 

least  from  the  true  value.  The  objective  analysis  estimated  value,  fgt,  at  point  g  and  time  t  com¬ 
puted  using  n  nearby  observations,  (/  =  1, 2, ....  n) ,  is  a  hnear  combination  of  the  preliminary 
field  estimate  and  the  deviations  of  the  observations  from  the  preliminary  field  estimate: 


j  =  1 


(9) 


The  linear  coefficients,  a.  (/  =  1, 2, n) ,  are  chosen  so  as  to  minimize  the  variance  of  the  field 
estimate  fgt: 

Vg,  =  <(4, -/*,)">  (10) 


=  <  (r;,)  ">  +  £  £  a, a.  ( </;,/' )  +  (e,,/;,)  +  (f'fy)  +  - 

1=  V=l 

i=  1 

Under  the  assumption  that  the  observational  errors  and  the  deviations  f^'  are  uncorrelated,  the 
variance  reduces  to 


i=  iy=i 


(11) 


i=  1 
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For  purposes  of  brevity,  the  symbol  d^j  =  will  be  introduced,  so  Equation  (11)  can  be 

rewritten  as 


v„  =  m 


n  n  n 

ss*  L  E -^E"!”*'- 

i=lj=l  /=! 


(12) 


A  necessary  condition  for  minimizing  with  respect  to  the  weights  is  that  for  =  1, 2, ....  n 


3a, 


_ 


=  2 


Li  =  1 


=  0,  or 


(13) 


J^a.(m,.^  +  d.^)  = 


(14) 


/  =  1 


In  order  to  solve  the  linear  system  (14),  the  covariances,  and  must  be  known  for 

i,k  =  1,2,  The  computer  code  described  in  Section  4  assumes  that  the  covariances  are 
functions  of  the  form 


C 


(15) 


where  aj  >  0  ^nd  bj  <  0-  The  parameter  is  an  a  priori  estimate  of  the  variance  of  the  prelimi¬ 
nary  field  estimate,  while  the  parameters  02  and  depend  on  the  spatial  and  temporal  separation 
between  points  of  interest.  One  choice  of  parameters  is 


r^.7)  ^ 

Uo  ) 

<  ^0  y 

(16) 


and 


1 

2 


(17) 


The  quantities  Ax^j,  Ay.j,  and  At.j  are  the  x,  y  and  t  components,  respectively,  of  the  distance 
between  positions  corresponding  to  field  values  and  fj.  The  parameters  Xq  and  yg  determine 
the  zero  crossing  of  the  correlation  function.  The  parameters  y^,  and  are  the  decay  scales 
for  the  correlation  function. 
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The  variance  of  the  objective  analysis  field  estimate,  which  is,  in  fact,  the  minimum  variance,  may 
be  found  by  substituting  Equation  (14)  into  Equation  (12): 


n 


-  m 


gg 


- 1 


(18) 


Equations  (9),  (14)  and  (18)  can  be  conveniently  expressed  in  vector-matrix  notation,  following 
Reference  [b].  Let 


A/= 

Q,  —  **•*  ’ 


c,= 


//111  +  dll  ... 


•••  ^nn'^^nnj 


^ gf  ~  (^gl>  •••>  ^g«)  ’ 


Then  Equations  (9),  (14)  and  (18)  become 


(19) 

(20) 


(21) 


(22) 


fgt  =  fgt  +  (a-  ^  (23) 

a  =  Cgfj\snd  (24) 

=  Var  (/,,)= /n^^-C^;^cJ.  (25) 

Note  that  the  subscript  g  refers  to  the  particular  grid  point  at  which  the  field  is  being  estimated, 
while  the  subscript  /  refers  to  the  field  being  estimated. 


3.2.2  Vector  Objective  Analysis 

Objective  analysis  can  also  be  used  to  estimate  vector  fields  such  as  velocity,  for  which  the  differ¬ 
ent  components  of  the  field  are  cross-correlated. 

In  vector  objective  analysis,  the  definitions  given  by  expressions  (19)  -  (22)  must  be  replaced  by 
analogous  versions,  in  which  vector  and  matrix  entries  are  replaced  by  matrices.  Equations  (23), 
(24)  and  (25)  remain  essentially  the  same. 
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In  vector  objective  analysis,  the  definitions  given  by  expressions  (19)  -  (22)  must  be  replaced  by 
analogous  versions,  in  which  vector  and  matrix  entries  are  replaced  by  matrices.  Equations  (23), 
(24)  and  (25)  remain  essentially  the  same. 


For  purposes  of  illustration,  the  equations  for  a  2-component  field,  U 
given.  Let 


u 

V  V  , 


will  now  be 


U: 


^it 


=  the  true  vector  field  value  at  position  i  at  time  t , 


(26) 


TTofci  _ 


^Vf^bsj 


=  the  observed  vector  field  value  at  position  i  at  time  t,  (27) 


En 


\^itJ 


=  observational  error  vector  at  position  i  and  time  t,  and 


(28) 


Uf, 


vvKy 


=  the  preliminary  field  estimate,  where 


UP,  = 


E  <v/,^ 


W  =  1 


/I  hi 
;■=  1 


and 


vfr  = 


I  <v/.. 


(29) 


(30) 


The  { l^j}  in  Equations  (29)  and  (30)  are  weights  as  in  Equation  (6).  Using  the  above  definitions, 
the  vector  equations  corresponding  to  Equations  (19)  through  (22)  may  be  given  as: 

Dl/  =  (t/ff^  -  t/f,, ...,  -  UP,) ,  (31) 

A  =  (Aj,  ...,A„) ,  where  (32) 
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Cjj  - 


^nl+^nl 


M  +  D 

nn  nn 


,  where 


M..  =  f  and 


(h.Vl  and 


Cgu  —  (M^j, ...,  - 


Now  the  vector  equations  analogous  to  Equations  (23),  (24)  and  (25)  become 


U  =  *  =  the  objective  analysis  field  estimate, 

=  U^g,  +  A-AU, 


A  -  C  fjCu  , 


=  covariance  matrix  for  Ug, 

(  {{Ug-  -  uj  (v^,  -  v^,)  >  "| 


gt  '^gtf  y'gt  ''gtf 


<  (“.r  -  -  y,t)  >  <  -  ^/)  > 


\\\'*gt  ‘*gt' ^’’gt  'gt 

^  ^gg~^gu^u^Ju- 


Equations  (26)  through  (39)  specify  the  vector  objective  analysis  estimates  for  the  field  compo¬ 
nents  and  their  covariances. 
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For  the  purposes  of  computer  coding,  it  is  more  convenient  to  re-express  Equations  (37)  and  (38)  so 
that  the  field  components  u  and  v  are  separated.  Given  the  definitions 


then  it  follows  that: 


M, 


m.CD 

V 

V 

m.^) 

y 

y 

di.) 

,il3, 

di4) 

Igt 


V  Vgr  J 


p 

Ugt 


V  v'gr  J 


Am, 


Am. 


Av, 


V  ”  y 


,  where 


^21  ■■■^2n 


^31- ■■^2n  ^41- ••^4n 


--1 


CgU  = 


CguCu 

,  and 

"’gP  • 

^  ^gP- 

•'"i? 

mjf) . 

Cu  = 


mlP  +dpp...mp'>  +dlP 


mpp+dpp...m(p+dpp 


mIP  +d[fK..m(p  +dlP 


mIP  +dip...mip  +dlP 

m^P  +dlp...m^P  +dlP 

+dltV-^nn  +^hn 


and 


^gt  =  M^g-CguCuCju- 
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4.0  COMPUTER  IMPLEMENTATION  OF  OA  ALGORITHM 


In  this  section,  we  describe  a  computer  implementation  in  FORTRAN  of  the  vector  OA  algo¬ 
rithms  presented  in  Section  3. 

A  scalar  version  of  the  OA  algorithm  in  FORTRAN  is  contained  in  the  Harvard  University 
GULFCAST  system.  The  code  was  modified  in  this  project  to  accept  vector  components. 

Program  Input.  The  velocity  observation  data  is  shown  in  Table  4-1.  The  OA  program  is  con¬ 
trolled  by  the  parameters  which  are  shown  in  Table  4-2. 


Table  4-1.  Velocity  Observation  Data 


Program 

Parameter 

Description 

DEPTH 

Depth  at  which  the  observation  was  made. 

TIM 

Time  at  which  the  observation  was  made. 

LON 

Longitude  at  which  the  observation  was  made. 

LAT 

Latitude  at  which  the  observation  was  made. 

U 

Zonal  component  of  velocity. 

V 

Meridional  component  of  velocity. 

UVAR 

Variance  of  the  zonal  component  of  velocity. 

WAR 

Variance  of  the  meridional  component  of  velocity. 

Program  Flow.  The  vector  objective  analysis  program  begins  by  reading  an  input  file  containing 
the  parameters  listed  in  the  Table  4-2.  The  last  line  of  the  parameter  file  contains  the  name  of  the 
file  containing  the  observation  data. 

The  observation  data  file  is  read,  and  the  data  is  stored  in  arrays.  Next  the  positions  of  the  grid 
points  in  the  interpolation  region  are  calculated.  The  following  steps  are  then  repeated  for  each 
grid  point. 


11 


Applied  Mathematics,  Inc. 


Table  4-2.  OA  Parameters 


Program 

Parameter 

Description 

IM 

Number  of  grid  points  in  the  zonal  directioa 

JM 

Number  of  grid  points  in  the  meridional  direction. 

DX 

Zonal  grid  spacing. 

DY 

Meridional  grid  spacing. 

CLNG 

Domain  centroid  longitude  (west  values  are  negative). 

CLAT 

Domain  centroid  latitude  (south  values  are  negative). 

XOFF 

Zonal  offset  for  the  domain  centroid. 

YOFF 

Meridional  offset  for  the  domain  centroid. 

THETAD 

Domain  rotation  (counterclockwise  from  east). 

NNCE 

Maximum  number  of  influential  observations  used  in  estimating  the  field  at  a  point  of 
interest. 

RNCE 

Spatial  radius  of  influence  of  observations. 

TNCE 

Temporal  radius  of  influence. 

XZERO 

2^ro  crossing  m  the  zonal  direction  of  the  analytical  correlation  function  (corresponding  to 
Xq  in  Equation  (16)). 

YZERO 

Zero  crossing  in  the  meridional  direction  of  the  analytical  correlation  function  (corre¬ 
sponding  to  yQ  m  Equation  (16)). 

XDCAY 

Zonal  decorrelation  scale  of  the  analytical  correlation  function  (corresponding  to  in 
Equation  (17)). 

YDCAY 

Meridional  decorrelation  scale  of  the  analytical  correlation  function  (corresponding  to 
in  Equation  (17)). 

TDCAY 

Temporal  decorrelation  scale  of  the  analytical  correlation  function  (corresponding  to  in 
Equation  (17)). 

UPHSE 

Zonal  phase  speed. 

VPHSE 

Meridional  phase  speed. 

OADAY 

Time  for  which  the  field  is  estimated. 
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A  subroutine,  SELECT,  searches  the  arrays  of  input  data  to  identify  observations  within  the  influ¬ 
ential  spatial  radius  (parameter  RNCE)  and  within  the  influential  temporal  radius  (parameter 
TNCE)  of  the  grid  point  of  interest.  If  more  than  the  user  specified  maximum  number  of  influen¬ 
tial  observations  (parameter  NNCE)  are  found,  the  NNCE  observations  that  are  most  highly  cor¬ 
related  with  the  point  of  interest  are  used. 

A  subroutine,  REMAV,  calculates  the  mean  and  variance  of  the  velocity  components  of  the  influ¬ 
ential  points  as  determined  by  subroutine  SELECT.  Subroutine  REMAV  also  subtracts  the  mean 
from  the  velocity  components,  thus  producing  the  vector  of  deviations  from  the  mean  which  will 
be  weighted  to  produce  the  estimated  velocity  components  at  the  point  of  interest. 

The  objective  analysis  calculations  are  performed  in  subroutine  OB  JAN.  This  subroutine  solves  a 
linear  system  for  the  weights  to  be  applied  to  the  deviations  calculated  by  REMAV.  Subroutine 
OBJAN  calls  subroutine  SETA  which  calculates  the  correlation  matrix  of  observations.  The  cor¬ 
relation  matrix  is  normalized  by  the  variances  of  the  velocity  components  of  the  influential  obser¬ 
vation.  This  matrix  is  block  diagonal  and  symmetric.  The  upper  left  block  represents  the 
correlation  between  the  first  components  of  velocity  while  the  lower  right  block  represents  the 
correlation  between  the  second  components  of  velocity.  The  remaining  off  diagonal  blocks  are 
zero  because  of  the  assumption  that  the  zonal  and  meridional  components  of  velocity  are  uncorre¬ 
lated.  Before  returning  control  to  OBJAN,  SETA  inverts  the  correlation  matrix  and  checks  for  ill- 
conditioning.  Subroutine  OBJAN  multiplies  the  inverted  correlation  matrix  by  the  vector  of  cor¬ 
relations  between  the  point  of  interest  and  the  observation  points.  The  resulting  vector  of  weights 
is  multiplied  elementwise  by  the  vector  of  deviations  from  the  mean  calculated  previously  by 
REMAV.  Once  the  previously  removed  velocity  mean  is  added  back,  the  estimation  of  the  veloc¬ 
ity  components  at  the  point  of  interest  is  complete.  The  algorithm  then  proceeds  to  the  calculation 
of  the  velocity  field  at  the  next  grid  point. 

Program  Output.  The  OA  program  is  run  for  the  specified  solution  time.  At  each  grid  point,  a 
velocity  field  and  associated  uncertainties  are  calculated.  The  uncertainty  in  each  velocity  compo¬ 
nent  at  a  grid  point  is  calculated  as  a  relative  quantity  with  respect  to  the  sample  variance  of  the 
input  values  determined  by  the  subroutine  SELECT.  These  relative  variances  are  converted  to 
absolute  error  values.  The  velocity  field  is  smoothed  using  a  Shapiro  filter  as  described  in  Refer¬ 
ences  [c]  and  [d].  This  filtering  smooths  regions  of  discontinuity  which  arise  as  the  influential 
points  change  during  the  OA  analysis.  The  latitude,  longitude,  x  component  of  velocity  and  vari¬ 
ance,  and  y  component  of  velocity  and  variance  are  stored  for  each  of  the  specified  grid  points. 
These  values  determine  the  velocity  and  error  fields  shown,  for  example,  in  Figure  6-3. 
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5.0  DATASET 


The  United  States  Coast  Guard  Research  and  Development  Center  provided  a  dataset  of  drifting 
buoy  positions  versus  time  (Reference  [e])  for  use  in  this  project.  Arthur  Allen  of  the  USCG  R&D 
Center  directed  the  deployment  of  several  different  types  of  drifter  buoys  from  the  Canadian 
Coast  Guard  Ship  Sir  Humphrey  Gilbert.  The  test  was  conducted  in  November  and  December, 
1993,  off  the  east  coast  of  Ae  Avalon  Peninsula  near  St.  John’s,  Newfoundland.  The  test  area, 
approximately  100  nm  on  a  side,  is  shown  in  Figure  5-1. 


5.1  DESCRIPTION 

A  MESflMET  buoy  was  moored  during  the  test.  This  buoy  was  equipped  with  sensors  for  air  and 
water  temperature,  wind  speed  and  direction,  air  pressure,  wave  height  and  period,  and  current. 
The  meteorological  instruments  were  mounted  on  a  mast  three  meters  above  the  sea  surface.  The 
MINIMET  buoy  was  anchored  in  water  approximately  600  feet  deep  on  a  1900  foot  tine  (which 
allowed  the  buoy  to  drift  in  a  circle  around  the  anchor  point).  Thirty  feet  below  the  buoy  was  a 
rectangular  window  shade  drogue  (2m  by  10m).  A  S4  current  meter  was  suspended  at  Im  depth 
from  a  float  70  feet  from  the  main  MINIMET  buoy.  The  MINIMET  buoy  served  as  a  reference 
point  for  the  deployment  of  the  drifter  buoys.  The  MINIMET  buoy  position  was  determined  every 
half  hour  using  GPS.  The  buoy  also  recorded  all  collected  data  on-board  in  EPROM  which  was 
read  when  the  buoy  was  recovered. 

In  addition  to  the  MINIMET  buoy,  ten  freely  drifting  buoys  were  deployed  at  various  locations 
and  times  during  the  test.  The  positions  of  these  buoys  were  collected  and  reported  by  the  Argos 
satellite  system.  The  initial  position  of  each  drifter  buoy  was  taken  to  be  the  deploying  ship’s  GPS 
position.  After  deployment,  the  drifter  buoys’  positions  were  calculated  and  recorded  by  the 
Argos  system.  The  Argos  system  uses  the  Doppler  shift  of  radio  waves  broadcast  at  frequency 
401.65  MHz  to  measure  the  position  of  transmitters  in  the  drifting  buoys. 

Service  Argos  provides  a  discrete  measure  of  positional  accuracy  called  the  “fix  quality”.  The 
qualities  are  given  in  Table  5-1. 


Table  5-1.  Argos  Position  Error 


Fix  Quality 

Positional  Error 

3 

<150m 

2 

<  350m 

1 

<  1000m 

0 

>  1000m 
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The  errors  given  in  Table  5- 1  are  one  standard  deviation  values  For  example,  68%  of  positions 
reported  with  fix  quality  3  are  assumed  to  be  within  150m  of  the  true  position.  The  GPS  positions 
are  considered  accurate  to  within  100m.  Any  GPS  positions  recorded  in  the  data  were  indicated 
with  a  fix  quahty  of  9. 

Figure  5-2  shows  the  time  periods  in  which  the  buoys  were  deployed.  Table  5-2  lists  the  times, 
buoy  ID  number  and  fix  qualities  for  each  reported  position.  Figure  5-3  presents  a  histogram  of 
the  buoy  fix  quahties. 

Three  different  types  of  drifter  buoys  were  deployed:  three  self  locating  datum  marker  buoys 
(SL-DMB)  manufactured  by  Seimac  Limited  (ID  numbers  2619, 4521, 4529),  five  Davis/CODE 
drifters  (ID  #4531  and  4532  manufactured  by  Clearwater  Instruments  and  ED  #4533, 4534, 4535 
manufactured  by  Technoceans),  and  two  WOCE  drifter  buoys  (ID  #  2590  and  2591  manufactured 
by  Clearwater  Instruments). 

Three  of  the  drifting  buoys  (2591, 4521, 4529)  were  deployed  at  the  location  of  the  MINIMET 
buoy.  Because  the  MINIMET  moved  several  miles  south  before  the  end  of  the  test,  not  all  of  these 
buoys  were  deployed  at  the  same  location. 

The  drift  paths  are  shown  in  Figure  5-4  through  5-6  separately  for  each  type  of  buoy  deployed. 
The  green  asterisk  denotes  the  buoy  initial  points;  the  red  asterisk  indicates  the  buoy  final  points. 
No  quality  zero  fix  positions  are  plotted. 


5.2  DETERMINATION  OF  BUOY  VELOCITY 


The  velocities  of  the  buoys  were  calculated  using  the  position,  time,  and  fix  quality  information 
provided  in  the  dataset.  The  buoys  were  assumed  to  move  on  a  constant  course  and  speed  track 
between  observation  points.  The  positional  data  were  filtered  to  remove  all  observations  with  fix 
quality  0,  and  the  great  circle  distances  and  bearings  between  successive  positions  were  calcu¬ 
lated.  The  average  velocities  of  the  buoys  were  determined  from  these  distances  and  elapsed 
times.  These  average  velocities  were  resolved  into  perpendicular  components  using  the  bearing 
angles.  The  position  and  time  of  the  first  observation  in  each  pair  were  used  as  the  position  and 
time  of  each  derived  velocity. 
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28  NOV  93 

Buoy  ID  I _ ^ 

Number  I 


1  DEC  93 


6  DEC  93  U  DEC  93 

I - 1 - 1 - 1 - 1 - 1 - 1 - 1- 


16  DEC  93 


4535 

4534 

4533 

4532 

4531 

4529 

4521 

2619 

2591 

2590 


- 1 

H 


Figure  5-2.  Buoy  Deployment  Time  Line 
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Table  5-2.  Buoy  Position  Hme,  ID  Number  and  Fix  Quality 


Buoy  ID  Fix 
#  Quality 


DTG 

Buoy  ID 
# 

Fix 

Quality 

i 

DEC  93 

041030Z 


1831Z 

2009Z 

2011Z 

2015Z 

2151Z 

2331Z 

2332Z 

050647Z 


0648Z 


1006Z 


lOlOZ 


lOlOZ 


1145Z 


1326Z 


1817Z 

1817Z 

958Z 

959Z 

il30Z 

il31Z 

2309Z 

[312Z 


0957Z 


lOOOZ 


DTG 

Buoy  ED 

Fix 

Quality 

J 

DEC  93 

061807Z 

2619 

1 

1807Z 

4521 

1 

1945Z 

2 

1946Z 

2619 

1 

2108Z 

1 

2108Z 

1 

2249Z 

2 

2249Z 


070623Z 


0805Z 


0947Z 


1933Z 


933Z 


2048Z 

2048Z 

2229Z 

080752Z 

0933Z 

0933Z 

1040Z 


1040Z 


221Z 

221Z 


090739Z 


0740Z 


**  0920Z 


0924Z 


1020Z 
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Tkble  5-2.  Buoy  Position  Time,  ID  Number  and  Fix  Quality  (Continued) 


Buoy  ID  Fix 
#  Quality 


i 

7EC93 

091202Z 

2619 

0 

1202Z 

4521 

1 

1339Z 

2619 

0 

1616Z 

4533 

9 

1647Z 

4531 

9 

1711Z 

4529 

9 

1731Z 

4521 

0 

1732Z 

4534 

9 

1804Z 

4535 

9 

1840Z 

4532 

9 

*  2005Z 

4532 

2 

*  2006Z 

4531 

2 

2049Z 

2619 

0 

2049Z 

4534 

2 

205 IZ 

4529 

0 

205 IZ 

4531 

2 

2051Z 

2052Z 

4521 

0 

2052Z 

4533 

2 

4521 

2 

2148Z 


2324Z 

2324Z 

2325Z 

2325Z 

2325Z 

2326Z 

■ 

■ 

0727Z 

0729Z 


0908Z 


Buoy  ID 


Fix 

Quality 


Buoy  ID  Fix 
#  Quality 


DEC  93 


2619 

0 

4529 

2 

4531 

2 

4532 

3 

4533 

2 

4534 

2 

4535 

3 

4531 

2 
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Table  5-2.  Buoy  Position  Time,  ID  Number  and  Fix  Quality  (Continued) 


DIG 

Buoy  ID 
# 

Fix 

Quality 

)£C95 

141154Z 

4531 

2 

1154Z 

2 

1155Z 

2 

1330Z 

4531 

2 

1339Z 

2590 

9 

1635Z 

4534 

9 

1702Z 

4529 

9 

Bggil 

9 

2 

1834Z 

9 

1913Z 

9 

150947Z 

0 

1900Z 

9 

HOI 

0 

202243Z 

IHQI 

240937Z 

2619 

0 

1932Z 

2619 

0 

DTG 

Buoy  ID 

Fix 

Quality 

DEC  93 

122041Z 

4531 

3 

2041Z 

4532 

3 

2041Z 

4533 

2 

130000Z 

4532 

2 

OOOlZ 

4529 

0 

OOOlZ 

4531 

2 

0653Z 

4532 

3 

1012Z 

4531 

2 

2590 

1 

4529 

■■ 

1033Z 

4532 

3 

1033Z 

■sail 

2 

1212Z 

4533 

2 

1212Z 

4534 

3 

1214Z 

1 

1214Z 

2 

mEm 

4531 

2 

4531 

2 

1821Z 

2590 

1 

1823Z 

4531 

3 

2000Z 

4531 

2 

2002Z 

2590 

2 

2157Z 

4531 

2 

2157Z 

4535 

2 

2336Z 

4535 

2 

2338Z 

4534 

2 

2339Z 

4531 

2 

2339Z 

4532 

2 

140638Z 

4529 

0 

0638Z 

4535 

2 

0640Z 

4533 

3 

lOOOZ 

4531 

3 

lOOOZ 

4535 

2 

lOOlZ 

4529 

0 

lOlOZ 

4529 

3 

lOlOZ 

4532 

3 

lOlOZ 

4533 

3 

lOlOZ 

4534 

3 

1013Z 

2590 

1 

1151Z 

4535 

2 

1154Z 

4529 

1 

DTG 

Buoy  ID 

Fix 

Quality 

L 

me  93 

112239Z 

4531 

2 

2239Z 

4533 

2 

2239Z 

4534 

2 

2240Z 

4532 

2 

2243Z 

4535 

2 

120015Z 

2590 

9 

■ESS3 

MEM 

2 

0 

3 

0703Z 

4533 

2 

0704Z 

2590 

1 

0704Z 

4531 

2 

0847Z 

2590 

2 

0847Z 

2 

0848Z 

2 

0848Z 

4532 

2 

*  1023Z 

4532 

2 

*  1023Z 

4533 

2 

1024Z 

4529 

0 

*  1024Z 

4531 

2 

1052Z 

4529 

0 

1053Z 

2619 

0 

1055Z 

4531 

2 

1055Z 

4532 

2 

1055Z 

4534 

2 

1058Z 

4533 

2 

1234Z 

4531 

2 

1234Z 

4533 

2 

1235Z 

4534 

2 

1237Z 

4529 

3 

1237Z 

4532 

2 

1653Z 

4529 

0 

1834Z 

2590 

2 

1834Z 

4531 

2 

1834Z 

4532 

2 

1834Z 

4533 

2 

1834Z 

4534 

2 

2015Z 

4529 

0 

2015Z 

4534 

2 

2041Z 

2590 

1 

2041Z 

4529 

0 
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Number  of  Fixes 


Fix  Quality 


Note;  Total  =  382 


Figure  5-3.  Histogram  of  Buoy  Fix  Quality 
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Latitude 


Longitude 


Figure  5-4.  SL-DMB  Buoy  Tracks 
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Longitude 


Figure  5-5.  Davis/CODE  Buoy  Tracks 
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Latitude 


Figure  5-6.  WOCE  Buoy  Tracks 
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5.3  VELOCITY  ERROR  MODEL 


The  velocity  variance  used  as  an  input  parameter  in  the  objective  analysis  was  estimated  using  the 
Argos  fix  position  and  the  associated  fix  quahty. 

An  Argos  fix  determines  a  buoy’s  position  at  a  specified  time  with  a  standard  deviation  given  by 
the  fix  quality.  Two  successive  Argos  fixes  of  a  buoy  determine  a  mean  velocity  with  error  as  fol¬ 
lows. 

Assume  that  two  successive  Argos  fixes  are  defined  by 


Each  fix  has  an  associated  fix  quality  a.,  i  =  1, 2.  The  uncertainty  in  position  y^  at  time  is 
assumed  to  be  a  circular  normal  distribution  with  standard  deviation  a.  =  ‘  =  An 

example  of  a  circular  normal  probability  distribution  with  a  =  1  is  shown  in  Figure  5-7 . 

The  buoy’s  mean  velocity  from  time  to  ^2  is 


.  \ 
X 

y. 


h  h 


Xo  -X, 


We  assume  that  the  Argos  fixes  are  independent.  The  uncertainties  of  the  components  of  velocity 
are  then  given  by 


a.  = 

X 


o.  = 

y 


h  h 


2  ,  2 
^1  +  ^2 
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Probability 


0.0  0.1  0.2  0.3  0.4 

Probability 


Figure  5-7.  Circular  Normal  Probability  Distribution  for  a  =  1 


26 


Applied  Mathematics,  Inc. 


6.0  OA  DERIVED  VELOCITY  FIELD 


The  objective  analysis  computer  code  discussed  in  Section  4  was  used  to  generate  surface  current 
velocity  field  estimates  for  the  dataset  described  in  Section  5.  The  OA  inputs  used  and  the  results 
are  presented  in  Section  6.1 

A  numerical  subsampling  experiment  conducted  to  examine  the  number  of  measurements 
required  to  represent  the  flow  field  in  a  typical  SAR  domain  is  described  in  Section  6.2. 

6.1  OA  USING  COAST  GUARD  BUOY  DATA 

Input  Parameters.  The  OA  input  parameters  used  with  this  dataset  are  presented  in  Table  6-1. 


Table  6-1.  OA  Input  Values 


Program 

Parameter 

Value 

Description 

IM 

23 

Number  of  grid  points  in  the  zonal  direction 

JM 

23 

Number  of  grid  points  in  the  meridional  direction 

DX 

2.25  km 

Zonal  grid  spacing  (km) 

DY 

3.0  km 

Meridional  grid  spacing  (km) 

CLNG 

051  48 W 

Domain  center  longitude  (dd.dd) 

CLAT 

47  30N 

Domain  center  latitude  (dd.dd) 

XOFF 

0.0  km 

Zonal  offset  for  the  domain  centroid 

YOFF 

0.0  km 

Meridional  offset  for  the  domain  centroid 

THETAD 

0.0  degrees 

Domain  rotation  (counterclockwise  from  east) 

NNCE 

10 

Number  of  influential  points 

RNCE 

100.0  km 

Radius  of  influence  (km) 

TNCE 

1. 0-3.0  days 

Influential  time  window  (days) 

XZERO 

80.0  km 

Correlation,  zonal  zero  crossing  (km) 

YZERO 

80.0  km 

Correlation,  meridional  zero  crossing  (km) 

XDCAY 

40.0  km 

Zonal  decorrelation  (decay)  scale  (km) 

YDCAY 

40.0  km 

Meridional  decorrelation  (decay)  scale  (km) 

TDCAY 

3.0  days 

Temporal  decorrelation  (decay)  scale  (km) 

UPHSE 

0.0  m/s 

Zonal  phase  speed 

VPHSE 

0.0  m/s 

Meridional  phase  speed 

OADAY 

343.977083 

Solution  time;  09  DEC  93, 23:27 

344.961806 

10  DEC  93, 23:05 

345.946528 

11  DEC  93, 22:43 

346.861806 

12  DEC  93, 20:41 
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The  domain  centroid  position  was  chosen  to  be  near  the  center  of  the  ocean  region  where  most  of 
the  buoys  drifted.  The  OA  grid  dimensions  and  spacing  were  also  chosen  to  contain  the  region 
which  was  dense  with  the  input  data  from  the  drifting  buoys.  The  buoy  track  density  is  shown  in 
Figure  6-1.  Data  with  fix  quahty  0  are  excluded  fi-om  Figure  6-1.  The  grid  covers  an  area  of  3267 

km^.  The  OA  region  is  a  rectangle  with  center  at  47  30N,  05 1  48W.  The  northern  and  southern 
edges  are  parallel  to  the  equator;  the  offset  point  and  domain  rotation  angle  are  all  zero.  Choice  of 
the  maximum  number  of  influential  points  is  determined  by  balancing  the  amount  of  computation 
required  with  the  desired  quahty  of  output.  Carter  and  Robinson  (Reference  [f])  state  that  using 
ten  influential  points  is  usually  sufficient.  The  radius  of  influence  was  set  to  100  km  so  that  obser¬ 
vations  from  the  entire  input  dataset  would  be  eligible  to  contribute  to  the  velocity  estimate  at 
each  grid  point.  This  maximizes  the  amount  of  data  available  for  use  in  computing  the  flow  field. 
The  influential  time  window  was  varied  between  one  and  three  days  to  demonstrate  the  effect  of 
accumulation  of  input  data. 

The  zero  crossing  and  decorrelation  scales  are  chosen  based  on  prior  knowledge  of  the  specific 
ocean  region.  They  are  related  to  a  quantity  known  as  the  Rossby  radius.  The  Rossby  radius  varies 
with  location  and,  in  general,  is  the  distance  at  which  the  earth’s  gravitational  forces  balance  with 
its  rotational  forces.  For  the  North  Atlantic  ocean  near  the  location  of  the  objective  analysis 
region,  the  Rossby  radius  is  40  km.  The  zero  crossing  distance  was  chosen  to  be  two  Rossby  radii. 
The  temporal  decay  scale  was  chosen  by  examining  the  current  meter  dataset  from  the  MINIMET 
buoy.  The  dominant  time  scale  of  current  activity  was  estimated  to  be  between  one  and  three  days; 
thus,  the  temporal  decay  scale  was  chosen  to  be  three  days.  Examining  the  derived  velocities  of 
the  drifting  buoys,  no  long  term  prevalent  motion  of  the  buoys  was  seen;  accordingly,  the  phase 
velocities  of  the  region  were  set  to  zero. 

Times  for  which  the  velocity  fields  were  computed  was  varied  from  Julian  days  343  to  346  (corre¬ 
sponding  to  09-12  DEC  1993)  to  correspond  with  the  time  period  in  which  buoys  were  deployed. 

Those  fixes  with  an  Argos  fix  quahty  equal  to  zero  were  not  used.  Additional  data  with  fix  quality 
greater  than  zero  was  omitted  if  there  were  multiple  fixes  of  a  buoy  obtained  during  a  single 
Argos  pass,  or  if  the  Argos  fix  was  inconsistent  with  the  surrounding  fixes.  Data  not  used  in  the 
OA  is  noted  in  Table  5-1,  with  indicating  a  multiple  fix  and  “**”  indicating  inconsistent  data. 

Figure  6-2  presents  a  graph  of  the  correlation  function  using  the  appropriate  values  from 
Table  6-1. 

OA  Output.  Time  dependent  OA  estimates  were  performed  for  the  four  day  period  09-12  DEC 
during  which  six  buoys  were  simultaneously  deployed  in  a  27.0  run  x  35.6  nm  region. 

The  nature  of  the  prevalent  flow  fluctuated  throughout  this  time  period.  We  accordingly  examined 
sensitivity  to  the  influential  time  window  by  performing  OA  estimates  with  an  influential  time 
window  of  one,  two,  and  three  days.  The  solution  time  used  for  each  estimate  was  the  latest 
ARGOS  fix  time  of  all  buoys  on  the  given  day.  All  data  from  time  090000Z  DEC  93  to  the  solu¬ 
tion  time  was  available  for  an  OA  estimate.  The  actual  amount  of  data  considered  was  determined 
by  the  influential  time  window.  The  solution  times  are  given  in  Table  6-2. 
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Loiituae 


Figure  6-1.  Buoy  Track  Density 
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Figure  6-2.  Correlation 
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Table  6-2.  OA  Solution  Times 


Solution  Time 

Day 

(hh:mm) 

09  DEC  93 

23:27 

10  DEC  93 

23:05 

11  DEC  93 

22:43 

12  DEC  93 

20:41 

Because  no  buoy  data  prior  to  09  DEC  93  was  used,  we  did  not  calculate  an  OA  estimate  with  the 
two  day  influential  time  window  on  day  09  DEC.  Similarly,  we  did  not  calculate  an  OA  estimate 
with  the  three  day  influential  time  window  on  the  days  09  DEC  and  10  DEC. 

Figures  6-3  through  6-20  are  plots  of  the  OA  velocity  field  estimates.  Each  estimate  is  displayed 
in  two  plots.  The  first  type  of  plot  (e.g.  Figure  6-3)  is  the  velocity  field  overlaid  on  the  absolute 
velocity  eiror  field  which  is  shown  in  color.  This  plot  includes  the  last  position  of  each  buoy,  plot¬ 
ted  as  red  asterisks,  prior  to  the  solution  time. 

The  second  type  of  plot  (e.g.  Figure  6-4)  presents  the  same  velocity  field  without  the  error  field 
and  the  last  velocity  estimate  for  each  buoy  prior  to  the  solution  time.  These  velocities  are  plotted 
as  either  red  or  blue  vectors.  The  vector  is  red  if  the  velocity  observation  is  within  four  hours  of 
the  solution  time,  otherwise  it  is  blue.  The  MINIMET  S4  current  data  was  available  for  the  days 
09  DEC  93  and  10  DEC  93.  The  MINIMET  current  estimate  is  plotted  as  a  green  vector  for  those 
days  when  it  was  available.  This  estimate  is  an  average  of  the  MINIMET  current  readings  imme¬ 
diately  prior  to  and  after  the  solution  time. 

A  buoy  velocity  estimate  requires  using  two  Argos  fixes.  The  buoy  data  plotted  in  Figures  6-3 
through  6-20  are  based  on  the  last  two  fixes  before  the  solution  time.  The  distance  between  these 
fix  positions  and  the  time  between  them  determine  a  mean  velocity  during  this  period.  The  first 
plot,  described  above,  displays  the  buoy  position  at  the  last  fix.  This  is  the  last  geographical  posi¬ 
tion  of  the  buoy  prior  to  the  solution  time.  The  second  plot  displays  the  buoy  velocity  at  the  next- 
to-last  fix  and  is  the  last  velocity  used  in  the  OA. 

The  buoy  positions  prior  to  the  solution  time  plotted  in  the  figures  indicate  where  the  latest  data  is 
located  that  was  used  for  the  OA  estimate.  The  OA  estimates  show  lower  uncertainty  in  the  vicin¬ 
ity  of  this  data  and  higher  uncertainty  in  regions  away  from  the  data.  The  influence  of  an  observa¬ 
tion  decays  with  distance  and  larger  uncertainty  occurs  as  estimates  are  extrapolated  from  the 
observed  data. 

The  last  velocity  estimate  of  each  buoy  prior  to  the  solution  time  plotted  in  the  figures  indicates 
the  latest  velocity  which  was  used  in  the  OA  estimate.  The  influence  of  an  observation  decays  in 
time,  and  the  older  data  plotted  in  blue  will  have  less  influence  than  those  plotted  in  red.  Some 
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Figure  6-3.  OA  Velocity  Field  and  Error  Field  Using  a  1  Day  Influential  Time  Window 
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Figure  6-4.  OA  Velocity  Field  Using  a  1  Day  Influential  Time  Window 
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Figure  6-5.  OA  Velocity  Field  and  Error  Field  Using  a  1  Day  Influential  Time  Window 
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Figure  6-6.  OA  Velocity  Field  Using  a  1  Day  Influential  Time  Window 
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Figure  6-7.  OA  Velocity  Field  and  Error  Field  Using  a  1  Day  Influential  Time  Window 
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Figure  6-8.  OA  Velocity  Field  Using  a  1  Day  Influential  Time  Window 
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Figure  6-9.  OA  Velocity  Field  and  Error  Field  Using  a  1  Day  Influential  Time  Window 
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Figure  6-10.  OA  Velocity  Field  Using  a  1  Day  Influential  Time  Window 
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Figure  6-11.  OA  Velocity  Field  and  Error  Field  Using  a  2  Day  Influential  Time  Window 
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Figure  6-12.  OA  Velocity  Field  Using  a  2  Day  Influential  Time  Window 
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Figure  6-13.  OA  Velocity  Field  and  Error  Field  Using  a  2  Day  Influential  Time  Window 
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Figure  6-14.  OA  Velocity  Field  Using  a  2  Day  Influential  Time  Window 
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Figure  6-15.  OA  Velocity  Field  and  Error  Field  Using  a  2  Day  Influential  Time  Window 
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Figure  6-16.  OA  Velocity  Field  Using  a  2  Day  Influential  Time  Window 
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Figure  6-17.  OA  Velocity  Field  and  Error  Field  Using  a  3  Day  Influential  Time  Window 
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Figure  6-18.  OA  Velocity  Field  Using  a  3  Day  Influential  Time  Window 
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Figure  6-19.  OA  Velocity  Field  and  Error  Field  Using  a  3  Day  Influential  Time  Window 
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Figure  6-20.  OA  Velocity  Field  Using  a  3  Day  Influential  Time  Window 
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variations  between  the  latest  buoy  velocity  prior  to  the  solution  time  and  the  estimated  velocity 
field  should  be  expected. 

The  flow  in  this  region  changes  significantly  from  day  to  day.  The  predominate  flow  in  the  region 
on  09  DEC  is  eastward  and  changes  direction  to  southwestward  on  10  DEC.  This  flow  then 
changes  to  northeastward  on  11  DEC,  and  is  northward  on  12  DEC.  The  one  day  influential  time 
window  plots  are  shown  in  Figures  6-3  through  6-10.  The  two  day  plots  are  shown  in  Figures  6-11 
through  6-16,  and  the  three  day  plots  are  shown  in  Figures  6-17  and  6-20.  The  sensitivity  to  the 
influential  time  window  can  be  seen  in  the  two  and  three  day  plots  on  days  10  DEC  and  11  DEC. 
These  velocity  estimates  contain  converging  and  diverging  regions  of  flow.  Such  regions  are  not 
part  of  the  flow  field  but  are  irregularities  due  to  the  inclusion  of  conflicting  data  in  the  OA  esti¬ 
mate.  However,  on  12  DEC  the  two  and  three  day  estimates  show  more  reahstic  fields  than  the 
one  day  estimate. 

This  analysis  demonstrates  that  the  choice  of  influential  time  window  is  an  important  parameter 
for  this  region  where  the  flow  field  is  evolving  on  scales  from  one  to  three  days.  It  is  important  to 
note  that  this  is  particularly  an  issue  when  data  is  extrapolated  from  the  region  containing  buoy 
data.  The  one,  two  and  three  day  plots  show  small  differences  in  the  vicinity  of  the  buoys.  How¬ 
ever,  as  data  is  extrapolated  from  this  area,  the  estimated  velocity  fields  exhibit  greater  differ¬ 
ences. 


6.2  SUBSAMPLING  EXPERIMENT 


A  subsampling  experiment,  using  objective  analysis,  was  conducted  in  order  to  examine  the  num¬ 
ber  of  measurements  required  to  represent  the  prevalent  flow  field  in  a  typical  S  AR  domain. 

A  test  case  domain  was  selected  which  is  110  km  by  110  km  and  centered  at  69.5°W,  38.5°N. 
With  a  resolution  of  5  km  this  domain  has  23  x  23  grid  points.  The  location  of  this  domain  is 
shown  in  Figure  6-21  with  respect  to  the  Harvard  Forecast  domain  (GULFCAST)  and  Coast 
Guard  District  1 . 

The  BioSynop  experiment  was  carried  out  during  Oct  -  Nov  1988  in  a  somewhat  larger  domain 
(600  X  300  km).  The  BioSynop  domain  is  marked  as  B  and  the  test  domain  is  marked  as  T  in  Fig¬ 
ure  6-21.  From  the  hydrographic  dataset  of  temperature  and  salinity,  the  dynamic  height  field  is 
obtained,  which  is  then  objectively  analyzed  for  the  test  domain.  The  geostrophic  velocity  com¬ 
ponents  are  then  obtained  as  the  spatial  gradients  of  the  streamfunction  field  (which  is  propor¬ 
tional  to  the  dynamic  height).  The  temperature  and  salinity  OAs  are  shown  in  Figures  6-22a  and 
6-22b.  The  dynamic  height  OA  is  shown  in  Figure  6-22c.  The  OA  input  values  are  shown  in 
Table  6-3.  The  u  and  v  components  are  shown  as  contour  plots  in  Figures  6-23a  and  6-23b.  The 
derived  vector  velocity  field  is  shown  in  Figure  6- 23c. 

The  major  oceanographic  feature  in  this  domain  is  a  cold  cyclonic  patch  of  water  in  the  mid-east¬ 
ern  region  of  the  domain.  This  domain  is  situated  just  north  of  the  Gulf  Stream,  in  the  surround¬ 
ing  recirculation  region,  generally  known  as  the  Slopewater  Gyre.  The  circulation  in  this  region 
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is  a  mesoscale  eddy  field  interacting  transportationwise  with  the  Gulf  Stream  through  sub-basin 
scale  interactions  (Reference  [g]). 

A  subsamphng  experiment  was  carried  out  using  this  dataset  in  the  test  domain  starting  from  the 
objectively  analyzed  field  shown  in  Figure  6-23c.  This  23  x  23  5  km  resolution  field  is  subsam- 
pl^  at  varying  resolutions  from  10  to  50  km.  The  resulting  fields  at  these  different  resolutions  are 
objectively  analyzed  on  the  original  5  km,  23  X  23  grid. 

The  vector  OA  fields  for  resolutions  10, 20, 35,  and  50  km  are  shown  in  Figures  6-24a,  6-25a, 
6-26a,  and  6-27 a,  respectively.  The  difference  vector  fields,  shown  in  Figures  6-24b,  6-25b, 
6-26b,  and  6-27b,  present  their  respective  differences  from  the  original  vector  field  (Figure  6- 
23c).  Note  that  the  difference  vectors  increase  with  increasing  grid  resolution. 

The  root  mean  square  (rms)  difference  error  vs.  resolution  is  shown  in  Figure  6-28.  The  rms  error 
is  seen  to  be  only  1  cm/sec  for  this  field  which  has  a  range  of  velocity  magnitude  between  ±30 
cm/sec  up  to  the  resolution  of  35  km,  beyond  which  the  rms  error  rises  to  5  -  10%  of  the  maxi¬ 
mum  velocity  magnitudes. 

These  subsampling  results  are  encouraging.  Even  with  just  sixteen  (4x4)  data  points  for  the  35 
km  resolution  in  the  S  AR  domain  (110x110  km),  the  OA  field  (Figure  6-26a)  depicts  most  of 
the  oceanographic  features  as  seen  by  very  high  resolution  samphng  (Figure  6-23c).  This  analysis 
suggest  a  buoy  deployment  of  (4  x  4)  or  even  a  strategically  placed  (3x3)  array  for  a  1 10  x  1 10 
km  SAR  region  with  a  35  km  resolution  would  suffice  to  map  this  region  using  objective  analysis. 


Table  6-3.  OA  Input  Values  for  Subsampling  Experiment 


Values 

Description 

69.5  W 

Domain  center  longitude  (dd.dd) 

38.5  E 

Domain  center  latitude  (dd.dd) 

5 

Zonal  grid  spacing  (km) 

5 

Meridional  grid  spacing  (km) 

80.0 

Correlation,  zonal  zero  crossing  (km) 

80.0 

Correlation,  meridional  zero  crossing  (km) 

40.0 

Zonal  decorrelation  (decay)  scale  (km) 

40.0 

Meridional  decorrelation  (decay)  scale  (km) 

7 

Temporal  decorrelation  (decay)  scale  (days) 

10 

Number  of  influential  points 

200 

Radius  of  influence  (km) 

7 

Influential  time  window  (days) 

266 

Solution  day 
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Latitude 
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40 
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75  73  71  69  67  65 
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B  BioSynop  domain 
T  Test  domain 


Figure  6-21.  BioSynop  and  Subsampling  Test  Domain 
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Figure  6-22b.  Salinity  OA  Map 
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Figure  6-23a.  U  Component  OA  Map 
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Figure  6-23b.  V  Component  OA  Map 
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Figure  6-23c.  Vector  OA  Field:  5  km  Resolution 


Applied  Mathematics,  Inc. 


0 .320E^-00 
MAXIMUM  VECTOR 


Figure  6-24a.  Vector  OA  Field:  10  km  Resolution 
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Figure  6- 25a.  Vector  OA  Field:  20  km  Resolution 
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Figure  6-26a.  Vector  OA  Field:  35  km  Resolution 
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Figure  6-27a.  Vector  OA  Field:  50  km  Resolution 
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Grid  Resolution  (km) 


Figure  6-28.  RMS  Error  vs.  Grid  Resolution 
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7.0  COMBINING  VELOCITY  FIELDS 


In  this  section,  we  describe  two  approaches  for  combining  or  melding  regularly  gridded  velocity 
data  output  from  two  or  more  ocean  circulation  models  in  order  to  obtain  a  single  optimally  esti¬ 
mated  velocity  field. 

In  order  to  demonstrate  the  two  methodologies,  we  begin  with  a  parent  domain  P  and  a  subdo¬ 
main  A  as  shown  in  Figure  7-1.  The  OA  is  used  as  the  full  field  of  the  meridional  velocity  com¬ 
ponent,  described  in  the  Subsampling  Experiment  in  Section  6.2.  Figure  7-2a  shows  the  OA  field 
and  Figure  7-2b  shows  the  associated  error  when  domain  P  is  sampled  with  5  km  resolution. 

Subdomain  A  is  sampled  at  5  km  resolution  yielding  a  ( 14  x  12)  dataset.  This  can  be  considered 
to  be  the  representative  dataset  for  a  fine-grid  regional  numerical  model  output.  These  data  points 
are  indicated  in  Figure  7-1  as  small  x’s. 

Similarly,  domain  P  is  subsampled  at  35  km  resolution  yielding  a  (4  x  4)  dataset.  This  data  set  is 
representative  of  a  coarse-grid  general  circulation  numerical  model  output.  These  data  points  are 
marked  in  Figure  7-1  as  asterisks. 

If  we  use  objective  analysis  on  the  individual  datasets  separately  to  obtain  information  on  the  par¬ 
ent  grid,  P,  for  a  resolution  of  5  km,  we  will  have  fields  with  large  error  values  where  there  is 
missing  data  for  the  whole  domain.  Figure  7-3a  shows  the  5  km  OA  on  the  parent  domain,  P , 
using  the  data  in  domain  A .  The  associated  error  field,  shown  in  Figure  7-3b,  contains  large  error 
values  in  the  region  containing  no  data  points.  Similarly,  Figures  7-4a  and  7-4b  show  the  5  km 
OA  on  the  parent  domain,  P,  using  the  coarse-grid  dataset  and  its  associated  error  field,  respec¬ 
tively. 

Our  goal  is  to  reconstruct  the  OA  field  for  the  whole  region  P  so  that  it  matches  as  closely  as  pos¬ 
sible  the  original  field  shown  in  Figure  7-2a. 

Method  1.  Objective  Analysis.  One  method  is  to  pool  all  the  observations  from  both  datasets, 
and  carry  out  an  objective  analysis  via  a  standard  weighting  procedure.  This  method  of  using 
both  the  observations  of  domain  A  (at  5  km  resolution)  and  of  domain  P  (at  35  km  resolution) 
results  in  a  field  where  jagged  “seams”  are  formed  along  the  boundary  between  the  two  regions. 
This  is  seen  in  the  resulting  OA  field  shown  in  Figure  7 -5a.  The  seams  reflect  the  irregular,  non¬ 
smooth  behavior  of  the  combined  field,  which  in  turn  is  due  to  the  inconsistency  between  the  fine 
gridded  and  coarse-gridded  observations  in  the  boundary  region.  The  development  of  seams  is 
consistent  with  the  high  error  regions  along  the  boundary  between  the  regions  shown  in  the  asso¬ 
ciated  error  map  (Figure  7-5b). 
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Figure  7-1.  Grid  Resolution 
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Figure  7-2a.  TYue  Field 
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Figure  7-2b.  Error  Field  Corresponding  to  Figure  7-2a  with  5  km  Sampling  Resolution 
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Figure  7-3a.  OA  for  P  Domain  Using  Observations  from  Domain  A  Only 
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Figure  7-3b.  Error  Field  {PIA),  Corresponding  to  Field  Shown  in  Figure  7-3a 
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Figure  7-4a.  P  Domain  OA  Using  Coarse  Grid  (35  km)  Observations 


74 


Applied  Mathematics,  Inc. 


Figure  7-4b.  Error  Field  Corresponding  to  OA  Field  Shown  in  Figure  7-4a 
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Figure  7-5a. 


OA  Field  Using  Data  From  Domains  A  and  P 
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Figure  7-5b.  Error  Field 


spending  to  Field  Shown  in  Figure  7-5a 
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Method  2.  Weighted  Comhination.  The  second  method  of  combining  the  model  outputs  is  to 
use  objective  analysis  to  separately  extend  each  model  estimate  to  the  entire  domain  P .  At  each 
point  of  the  extended  grid,  take  a  weighted  combination  of  the  two  resulting  field  values.  The 

weights  are  determined  by  the  variances:  that  is,  if  models  1  and  2  give  field  estimates  fn  and  fi2 
at  grid  point  i  with  attached  variances  0^2  ’  respectively,  then  the  melded  field  estimate  at 

grid  point  /  is 


+ 


(40) 


Under  the  assumption  that  fn  and  fi2  are  independent,  the  estimate  fi 


nas  variance 


2  2 

2  ,  2“' 

0.1  + a.2 


(41) 


The  melded  OA  estimate  given  by  Equation  (40)  is  shown  in  Figure  7-6a.  The  modified  error 
field  as  given  by  Equation  (41)  is  shown  in  Figure  7-6b.  The  field  shown  in  Figure  7-6a  agrees 
with  the  field  in  Figure  7-2a,  and  does  not  have  any  unphysical  characteristics  like  the  seam  along 
the  boundary  between  the  two  regions. 
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Figure  7-6a.  Weighted  Combination  Field 
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Figure  7-6b.  Error  Field  Estimate  Corresponding  to  Field  Shown  in  Figure  7-6a 
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8.0  AREA  OF  UNCERTAINTY  PROJECTION  USING  OA 
VELOCITY  FIELD 


The  velocity  field  and  its  associated  error  field  can  be  used  to  project  forward  in  time  an  area  of 
uncertainty  (AOU).  This  projected  AOU  can  be  used,  for  example,  to  help  determine  where  to 
conduct  a  S  AR  mission. 

In  this  section,  we  present  a  methodology  for  calculating  the  AOU  projection  based  on  an  OA 
generated  velocity  and  error  field. 

Goal.  Project  an  AOU  forward  in  time,  given  a  velocity  field  which  has  uncertainty  at  each  point 
in  space. 

Method.  Choose  (Xg,  yg)  by  Monte  Carlo  from  the  original  AOU.  Project  forward  to  time  T  to 
get  the  trajectory.  Repeat  the  process  N  times,  to  get  an  AOU  at  time  T. 

Algorithm.  The  time  projection  also  uses  Monte  Carlo,  because  of  the  uncertainty  in  velocities. 
We  make  the  following  definitions: 

At  =  time  step  used  in  projection; 

(Xj,  yj)  s  projected  position  at  time  jAt; 

(Uj,  Vj)  s  OA  velocity  estimate  at  position  {Xj,  yj) ; 

(a^ ,  )  =  OA  velocity  variances  at  position  {Xj,  yj) ;  and 

(Cy ,  c^)  =  correlations  between  velocities  at  positions  {Xj,  yj)  and  j, 

{u  and  V  velocities  are  assumed  to  have  no  cross-correlation). 

According  to  the  OA  model  described  in  Section  3.0,  and  have  the  form 
c  =  (l-a2)e^\ 

where 
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e  ,  e  ,  6  ,  6  =  random  variables  to  be  defined  below. 


U-  V-  U-  V- 

J  J  J  J 


We  use  the  following  iterative  procedure  for  determining  a  Monte  Carlo  trajectory; 

where  £  ~  iV(0,  a  ) ;  e  ~  MO,  a  ) . 

“0  “o  ^0 

(^2’ >’2)  =  +At(Ui,Vi)  +  At(e^,e^), 

where 


£  =  c,  £  +6  ;6  ~M0,  \o^  -cfaf),and 

“0  1^0  “1  “1  ^  V  “)  “0  “0 


In  general, 


where 


£  =  c,  £  +6  ;  6  ~A^(0,  lo^  —c^  ), 

Uj  Uj-I  UJ.^  U/  Uj  ''  ’a/  Uj  Uj_i  Uj_/' 


and 


£=c„£  +6;6~  N{0,  lo^  -  ). 

V  V-  .  V  V  ’  V  '  ’  a/  V  V.  .  V  y 

j  V  - 1  j  - 1  j  J 


Note  that  £,  is  defined  so  that 


c^(e„)  =  ; 


<e,e.,  > 


“y  a.  a. 


S’ 


“j  “yi 


and  similar  equations  hold  for  £^ ;  thus,  velocity  deviation  variances  and  covariances  have  an 
appropriate  form. 
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Note  that  if  <  0,  we  may  simply  take  6^^  =  0;  and  similarly  for  6^ . 

Example.  Figure  8-1  shows  a  velocity  field  and  associated  error  field.  Figure  8-2  illustrates  a  12 
hour  time  projection  of  3  circular  AOUs  of  radius  2  nm  for  the  field  shown  in  Figure  8-1.  The 
green  region  denotes  the  initial  AOU;  the  red  region  the  projected  AOU.  A  total  of  288  points 
were  used  in  the  Monte  Carlo  for  this  display. 
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Figure  8-1.  Velocity  and  Error  Field  Used  for  AOU  Projection 
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Figure  8-2.  12  Hour  Time  Projection  of  3  Circular  AOUs  of  Radius  2  rnn 
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9.0  CONCLUSIONS 


In  this  section,  we  present  the  conclusions  of  the  Phase  I  effort. 

•  Objective  analysis  numerical  techniques  can  be  applied  to  Coast  Guard  observational 
data  to  provide  good  estimates  of  ocean  surface  velocity  field  and  associated  error. 

•  Using  a  TAC-3  (HP  700  series)  computer  system,  the  time  required  to  calculate  and  dis¬ 
play  ocean  surface  velocity  fields  for  typical  SAR  regions  is  less  than  10  seconds.  (Note: 
Code  has  not  been  optimized.  Run  time  will  be  considerably  reduced  with  optimized  code.) 

•  It  is  technically  feasible  to  develop  a  system  for  use  in  real  time  by  the  Coast  Guard  to 
provide  continuously  updated  ocean  surface  velocity  estimates  as  observed  data  are 
received. 
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